Response retention and apparent motion effect in visual cortex models

Apparent motion is a visual illusion in which stationary stimuli, flashing in distinct spatial locations at certain time intervals, are perceived as one stimulus moving between these locations. In the primary visual cortex, apparent-motion stimuli produce smooth spatio-temporal patterns of activity similar to those produced by continuously moving stimuli. An important prerequisite for producing such activity patterns is prolongation of responses to brief stimuli. Indeed, a brief stimulus can evoke in the visual cortex a long response, outlasting the stimulus by hundreds of milliseconds. Here we use firing-rate based models with simple ring structure, and biologically-detailed conductance-based refractory density (CBRD) model with retinotopic space representation to analyze the response retention and the origin of smooth profiles of activity in response to apparent-motion stimuli. We show that the strength of recurrent connectivity is the major factor that endorses neuronal networks with the ability for response retention. The same strengths of recurrent connections mediate the appearance of bump attractor in the ring models. Factors such as synaptic depression, NMDA receptor mediated currents, and conductances regulating spike adaptation influence response retention, but cannot substitute for the weakness of recurrent connections to reproduce response retention in models with weak connectivity. However, the weakness of lateral recurrent connections can be compensated by layering: in multi-layer models even with weaker connections the activity retains due to its feedforward propagation from layer to layer. Using CBRD model with retinotopic space representation we further show that smooth spatio-temporal profiles of activity in response to apparent-motion stimuli are produced in the models expressing response retention, but not in the models that fail to produce response retention. Together, these results demonstrate a link between response retention and the ability of neuronal networks to generate spatio-temporal patterns of activity, which are compatible with perception of apparent motion.


Model of V1.
V1 is modeled as a 2-dimensional continuum of neuronal populations.Each point of the cortical continuum contains two neuronal populations, excitatory (E) and inhibitory (I), connected by AMPA (α-amino-3-hydroxy-5-methyl-4-isoxazolepropionic acid), NMDA (N-methyl-d-aspartate), and GABA-A mediated synapses providing recurrent intracortical interactions and AMPA and NMDA for the geniculate input.The strengths of the external connections correspond to the pinwheel architecture, and thus neurons receive inputs in accordance with their orientation and direction preferences, taken from the experimental map shown in Fig. 5.The profile of the intracortical connections is isotropic, i.e., the maximum conductances depend on the distance between the preand postsynaptic populations.
We define a neuronal population as a group of similar neurons that receives both a common input from presynaptic populations and an individual noise.The mathematical description of each population is based on the probability density approach, namely, the CBRD approach [4], The model for each E-and I-neuronal populations takes into account two neuronal compartments and a set of voltage-gated ionic currents, including the adaptation currents.
The neuronal population firing rate determines the presynaptic firing rate, which in turn controls the dynamics of synaptic conductances.The presynaptic firing rate predetermined by the excitatory population firing rate determines the dynamics of AMPA and NMDA synaptic conductances.The inhibitory population controls the GABA conductance.The synaptic conductances are the input signals for the postsynaptic neuronal populations.The membrane voltage distribution across t  determines the output firing rate and so on.
Model of a single neuronal population.In the CBRD model, the membrane potentials and ionic channel states of the neurons of one population are dispersed due to the noise, and thus they are distributed in a space of neuronal refractority.The refractority state is characterized with the time elapsed since the last spike, t  .Single population dynamics are governed by the equations for neuronal density, the mean over noise realizations voltage, and the gating variables.Neurons that fire contribute to the population firing rate, which is the output measure of the population activity.
The firing rate can be quite precisely and computationally effectively calculated by solving a system of equations in partial derivatives, 1-d transport equations.The equations govern the evolution of neuronal states in the t  -phase space.They contain the Hodgkin-Huxley equations for the membrane voltage and gating variables, parameterized by  * , as well as the equation for the neuronal density in  * -space,   (,  * ).The output characteristic of the population's activity is the firing rate   (), which is equal to   in the state of a spike,  * = 0.The equations written below describe an excitatory population of adaptive regular spiking pyramidal cells according to [4] and [5].
In comparison with one-compartment model, the extra parameters is the ratio of dendritic to somatic conductances  and the dendritic length.We assume that the inhibitory synapses are located at soma, contributing into the somatic synaptic current   , whereas the excitatory synapses are at dendrites, determining the dendritic synaptic current   .Parameterized by  * , the governing equations for the population  ( or ) are as follows: ( ( where    (,  * ) is the total conductance;  is the square ratio of the dendritic length to the characteristic length.The somatic and dendritic synaptic currents   and   are calculated as where the differential operator represents the solution of the reverse problem of dendritic current estimation from somatically registered-like conductances [6].We imply that the synaptic conductance kinetics is estimated from somatic responses to stimulation of presynaptic neuronal population, thus it implicitly accounts not only the kinetics of synaptic channels but also the dendritic and axonal propagation delays.For the dendritic compartment, the differential operator sharpens the transient effect of the channels, thus providing better agreement between somatic postsynaptic currents and potentials.Note that we take into account this sharpening only for glutamatergic channels by placing them on the dendritic compartment [7].
Hazard function.The source term in the eq.( 1) is the hazard function  which is defined as the probability for a single neuron to generate a spike, if known actual neuron state variables.The hazard function  has been approximated in [5] for the case of color noise as a function of   () and   (,  * ),, and the noise amplitude in the resting state   0 , the spike threshold voltage  ℎ and the ratio of membrane to noise time constants  =   /  : where  is the membrane potential relative to the threshold, scaled by noise amplitude;  is the hazard for a neuron to cross the threshold because of noise, derived analytically and approximated by exponential and polynomial for convenience;  is the hazard for a neuron to fire because of depolarization due to deterministic drive, i.e. the hazard due to drift in the voltage phase space.Note that the  -function is independent of the basic neuron model and does not contain any free parameters or functions for fitting to any particular case.Thus,  -function is the same for excitatory and inhibitory populations.
Voltage-dependent channels of excitatory neurons.The ionic currents  −  (  , ,  * ) = −  −   −   −   include the voltage-dependent potassium currents   and   responsible for spike repolarization, the slow potassium current   that contributes to spike frequency adaptation, the cation current   and the potassium current   , implicitly dependent on calcium dynamics, which also contributes to spike frequency adaptation.The approximating formulas for the currents   ,   ,   ,   and   are taken from [8]; the approximation for   is from [9].
where  *  is such that  Here   0 is the total somatic conductance at rest, and   is the total synaptic conductance;  is the membrane area.The dependence of  ℎ ( * ) is taken from a full single neuron model [4], allowing to take into account the effect of sodium channel inactivation on the threshold dynamics.  is the noise amplitude meaning the dispersion of individual neuron's voltage fluctuations in a stationary state.Its scaling with   approximately reflects the fact of the synaptic noise increase with the increase of mean synaptic drive.
For interneurons: When calculating the dynamics of a neural population, the integration of eqs.(12-26) determines the evolution of the distribution of voltage   across  * .Then, the effect of crossing the threshold and the diffusion due to noise are taken into account by  -function, eq.( 19), substituted into the equation for neuronal density (12).The integral (29) results in the output firing rate   ().

Representative neurons
Representative neurons of each of the populations were modeled according with the basic single neuron model with the same synaptic inputs as for the populations.The activity of the representative neurons does not affect the network.The representative neuron of E-population, for instance, is described by the equations for the membrane voltage, eq.(13,14), where the sum of partial derivatives were substituted by the total derivative in time t, and the sodium current was explicitly present in the right-hand part of Eq. ( 13).The sodium current dependent on voltage U was approximated by the 4-state Markov model [8]: Lognormal distribution of synaptic weights within each population.In order to introduce realistic, lognormal distribution of synaptic weights within a population  ( or ), the CBRD-approach has been generalized [11].In this case, instead of equal total synaptic current, neurons receive lognormally distributed current.For the current scaled by its mean across the distribution, , the distribution is The membrane potential of neurons parameterized with ,    , can be found as where    ( * ) is the unperturbed potential defined for zero synaptic input.
The density of neurons parameterized by  and distributed in the phase space  * is denoted as    (,  * ).Calculation of    (,  * ) requires solving of a continuum of eqs.(1) (or eq.( 18)) for    instead of   with (   ,    /).The output firing rate is defined as In numerical simulations, we set the parameter of the lognormal distribution   = 0.5 and discretized the  -space by 10 intervals.
Synaptic connections.The types of synapses are denoted by the types of mediator and postsynaptic neurons as follows: (, ), (, ) and (, ) with the postsynaptic index  =E orI.
The kinetics of the synaptic conductances is described following [7]   They are defined by the relations for the presynaptic rates  , () dependent on the somatic rates   , where presynaptic and postsynaptic populations are indexed by i and j, correspondingly.In the case each sending signals from either transient or sustained thalamic cells.The width and elongation of the footprints were 0.3 deg.and 0.8 deg., respectively.The neighboring V1 neurons that belong to different orientational hypercolumns have footprints of similar shapes and prefer the same orientation but opposite directions of the stimulus movement.In the T-S cell based mechanism of direction selectivity, the input firing rate is  ℎ  (, , ) =  ̃ ̃−1 (, ,  ̃,  ̃)  (,, ̃, ̃)( ̃,  ̃, ),
with the time constant    = 3.5  ,    = 3.5  , 5  instead of   ,   and   .Orientation and direction selectivity is determined by the footprint of the LGNto-V1 projections as well as by properties of the LGN neurons.The thalamic input is determined as the firing rate  ℎ  (, , ), which is a convolution of the LGN neuronal activity of the transient and sustained cells,   (, , ) and (, , ), respectively, with the footprint function  −1 (, ,  ̃,  ̃, ), where (, ,  ̃,  ̃) is the index of the LGN neuronal population, which attributes  or  indexes, respectively:  ℎ  (, , ) =  ̃ ̃∑  −1 (, ,  ̃,  ̃, (, ,  ̃,  ̃))  (,, ̃, ̃)( ̃,  ̃, ) with a second order ordinary differential equation, where the input is the presynaptic firing rate, i.e. as follows  ′, () =  ′,  ′, () ,  , () =  ,  , () ,Here  , is the presynaptic firing rate determined by axons of the population  on the postsynaptic population .In neglect of spatial propagation and temporal delays the presynaptic firing rate is equivalent to the somatic firing rate, i.e.  , ≡   .The index  is the synapse type,  = ,  or ; the index  = ℎ means thalamic input for  = ′;  =  for  =  or ; and  =  for  = ;  , is the maximum conductance,   , and   , are the rise and decay time constants.We imply that the synaptic time constants are estimated from the somatic responses to the stimulation of a presynaptic neuronal population, thus these time constants characterize not only synaptic channel kinetics but the dendritic and axonal propagation delays as well.The time scale  , is chosen in the form of Eq. (43) in order to provide independence of the maximum of  , () on   , and   , , when  , () is evoked by a short pulse of  , ().
is the magnesium (Mg2 +) concentration in mM;  , () is the non-dimensional synaptic conductance which is approximated by the second order ordinary differential equation: